\documentclass[utf8]{ctexart}

\usepackage[a4paper,left=1.25in,right=1.25in,top=1in,bottom=1in]{geometry}
\usepackage{listings}
\usepackage{graphicx}
\usepackage{subfigure}
\usepackage{booktabs}
\usepackage{amsmath}
\usepackage{amsthm}
\usepackage{amsfonts}
\usepackage{float}
\usepackage{indentfirst}
\usepackage{gnuplot-lua-tikz}
\usepackage{tikz}
\usetikzlibrary{shapes,arrows}
\usetikzlibrary{shapes.geometric, arrows}
\usepackage{algorithm}
\usepackage{algorithmic}
\usepackage{newclude}
\usepackage[perpage]{footmisc}

\graphicspath{ {images/} }
\raggedbottom	% 令页面在垂直方向向顶部对齐
\renewcommand\qedsymbol{QED}
\renewcommand{\figurename}{Figure}
\renewcommand{\proofname}{Proof}
\newcommand{\sign}[1]{\mathrm{sgn}(#1)}
\everymath{\displaystyle}   % 行内公式采用行间公式格式排列

\title{Numerical Analysis Homework \#6}
\author{Yin Wenliang\qquad 3200101893}
\date{}

\begin{document}
\maketitle
\CTEXsetup[format={\Large\bfseries}]{section}

\section{Theoretical questions}
\begin{enumerate}
\item
  % 问题1
  \begin{proof}
    \begin{enumerate}
    \item %(a)
      We consider the interpolation polynomial satisfying the condition required. Then the table of divided difference is as follows.
      \begin{table}[H]
        \centering
        \begin{tabular}{c|cccc}
          -1 & $f(-1)$ & & &\\
          0 & $f(0)$ & $f(0)-f(-1)$ & &\\
          0 & $f(0)$ & $f'(0)$ & $f'(0)-f(0)+f(-1)$ & \\
          1 & $f(1)$ & $f(1)-f(0)$ & $f(1)-f(0)-f'(0)$ & $\frac{f(1)-2f'(0)-f(-1)}{2}$\\
        \end{tabular}
      \end{table}
      Hence
      \begin{align*}
        p_3(y;-1,0,0,1;t) &= f(-1) + (f(0)-f(-1))(t+1)+ (f'(0)-f(0)+f(-1))t(t+1)\\
        &+\frac{f(1)-2f'(0)-f(-1)}{2}t^2(t+1).
      \end{align*}
      Then we have
      \begin{align*}
        &\int_{-1}^1p_3(y;-1,0,0,1;t)dt\\
        & = \int_{-1}^1[f(0)+\frac{f(-1)+f(1)-2f(0)}{2}t^2 + \frac{f(1)-2f'(0)-f(-1)}{2}t^3]dt\\
        & = \frac{1}{3}[f(-1)+4f(0)+f(1)] = \int_{-1}^1y(t)dt - E^S(y),
      \end{align*}
      which completes the proof.

    \item %(b)
      By (a) we have
      \begin{align*}
        E^S(y) &= \int_{-1}^1[y(t)-p_3(y;-1,0,0,1;t)]dt\\
        &= \int_{-1}^1\frac{f^{(4)}(\xi)}{4!}(t+1)t^2(t-1)dt\\
        &= \frac{f^{(4)}(\eta)}{24}\int_{-1}^1(t^4-t^2)dt = -\frac{f^{(4)}(\eta)}{90},
      \end{align*}
      where the second step follows from Theorem 2.35, the third step follows from the integral mean value theorem(Theorem C.71) and $\xi,\eta \in (-1,1)$.
    \item %(c)
      Firstly, we make a change of variable. Set $x=\frac{b-a}{2}t+\frac{b+a}{2}$, and we have
      \begin{align*}
        \int_{a}^by(x)dx &= \int_{-1}^1y(\frac{b-a}{2}t+\frac{b+a}{2})\frac{b-a}{2}dt\\
        &= \frac{b-a}{2}\int_{-1}^1p_3(y';-1,0,0,1;t)dt+E^S\\
        &= \int_a^bp_3(y;a,\frac{a+b}{2},\frac{a+b}{2},b;x)dx + E^S\\
        &= \frac{b-a}{6}[f(a)+4f(\frac{a+b}{2})+f(b)]+E^S.
      \end{align*}
      Hence,
      \begin{align*}
        E^S &= \int_a^b[y(x)-p_3(y;a,\frac{a+b}{2},\frac{a+b}{2},b;x)]dx\\
        &= \int_a^b\frac{f^{(4)}(\xi)}{4!}(x-a)(x-\frac{a+b}{2})^2(x-b)dx\\
        &= -\frac{(b-a)^5}{2880}f^{(4)}(\xi),
      \end{align*}
      where $\xi \in(0,1)$.\par
      Next, we apply above conclusion by subintervals and derive composite Simpson's rule
      \begin{equation*}
        I_n^s(f) = \frac{h}{3}[f(x_0)+4f(x_1)+2f(x_2)+4f(x_3)+2f(x_4)+\cdots+4f(x_{n-1}+f(x_n))],
      \end{equation*}
      where $h=\frac{b-a}{n}$,$x_k=a+kh$ and $n\in 2\mathbb{N}^+$.\par
      Sum up the errors, and we have
      \begin{align*}
        E_n^S(f) &= -\sum_{k=1}^{n/2}\frac{(2h)^5}{2880}f^{(4)}(\xi_k)\\
        &= -\frac{b-a}{180}h^4(\frac{2}{n}\sum_{k=1}^{n/2}f^{(4)}(\xi_k))\\
        &= -\frac{b-a}{180}h^4f^{(4)}(\xi),
      \end{align*}
      where the last step follows from the intermediate value Theorem C.39 and the fact $f^{(4)}\in \mathcal{C}[a,b]$. Here $\xi \in (a,b)$.
      \end{enumerate}
  \end{proof}
  
\item
  % 问题2
  \textbf{Solution}\par
  \begin{enumerate}
  \item %(a)
    By Theorem 6.17, the required condition yields
    \begin{equation*}
      \vert E_n^T(f)\vert = \vert -\frac{b-a}{12}h^2f''(\xi)\vert = \frac{1}{12n^2}\vert f''(\xi)\vert \leq 0.5\times 10^{-6}.
    \end{equation*}
    Notice that $f"(x) = (4x^2-2)e^{-x^2} < 2e^{-1},\forall x\in(0,1)$, hence
    \begin{equation*}
      n\geq \sqrt{\frac{2e^{-1}}{12\times 0.5\times10^{-6}}} \approx 350.18
    \end{equation*}
    Therefore 351 subintervals are required.
  \item %(b)
    Similarly, by Theorem 6.19, we have
    \begin{equation*}
      \vert E_n^S(f) \vert = \vert -\frac{b-a}{180}h^4f^{(4)}(\xi)\vert = \frac{1}{180n^4}\vert f^{(4)}(\xi)\vert \leq 0.5\times 10^{-6}.
    \end{equation*}
    Notice that $f^{(4)}(x) = (16x^4-48x^2+12)e^{-x^2}<12,\forall x\in(0,1)$,hence
    \begin{equation*}
      n\geq \sqrt[4]{\frac{12}{180\times0.5\times10^{-6}}} \approx 19.11.
    \end{equation*}
    Therefore 20 subintervals are required.
    \end{enumerate}
  
\item
  % 问题3
  \textbf{Solution}\par
  \begin{enumerate}
  \item %(a)
    Since $\mathbb{P}_1 = span(1,t)$, the statement that requires to be proved is equivalent to $<1,\pi_2(t)>=0$ and $<t,\pi_2(t)>=0$, which yield
    \begin{align*}
      \int_0^{+\infty}(t^2+at+b)e^{-t}dt &= 2+a+b = 0,\\
      \int_0^{+\infty}t(t^2+at+b)e^{-t}dt &= 6+2a+b = 0,
    \end{align*}
    respectively. Hence $a = -4,b = 2$, and the orthogonal polynomial is
    \begin{equation*}
      \pi_2(t) = t^2 - 4t + 2.
    \end{equation*}    

  \item %(b)
    The roots of $\pi_2(t)$ are $t_1 = 2-\sqrt{2},t_2 = 2+\sqrt{2}$. To calculate $\omega_1$ and $\omega_2$, consider
    \begin{align*}
      \omega_1 + \omega_2 &= \int_0^{+\infty}e^{-t}dt = 1,\\
      t_1\omega_1 + t_2\omega_2 &= \int_0^{+\infty}te^{-t}dt = 1,
    \end{align*}
    which yield
    \begin{align*}
      \omega_1 &= \frac{1-t_2}{t_1-t_2} = \frac{2+\sqrt{2}}{4},\\
      \omega_2 &= \frac{t_1-1}{t_1 - t_2} = \frac{2-\sqrt{2}}{4}.
    \end{align*}
    The desired two-point Gauss-Laguerre quadrature formula is thus
    \begin{equation*}
      I_2^G(f) = \frac{2+\sqrt{2}}{4}f(2-\sqrt{2}) + \frac{2-\sqrt{2}}{4}f(2+\sqrt{2}).
    \end{equation*}
    Hence we have
    \begin{equation*}
      \int_0^{+\infty}f(t)e^{-t}dt = I_2^G(f) + E_2(f),
    \end{equation*}
    with
    \begin{equation*}
      E_2(f) = \frac{f^{(4)}(\tau)}{4!}\int_0^{+\infty}e^{-t}(t^2-4t+2)^2dt = \frac{f^{(4)}(\tau)}{6},
    \end{equation*}
    where $\tau\in(0,\infty)$, and the first equation follows from Theorem 6.35.

  \item %(c)
    With $f(t) = \frac{1}{1+t}$, we have
    \begin{equation*}
      \int_0^{+\infty}f(t)e^{-t}dt = \frac{2+\sqrt{2}}{4}f(2-\sqrt{2}) + \frac{2-\sqrt{2}}{4}f(2+\sqrt{2}) \approx 0.5714286,
    \end{equation*}
    combining with $f^{(4)}(t) = \frac{24}{(1+t)^5} < 24, \forall t\in (0,\infty)$, there exists
    \begin{equation*}
      E_2(f) =\frac{f^{(4)}(\tau)}{6} < 4.
    \end{equation*}
    By using the exact value $I=0.596347361\cdots$, we have
    \begin{equation*}
      E_2^{true}(f) = I-\frac{2+\sqrt{2}}{4}f(2-\sqrt{2}) - \frac{2-\sqrt{2}}{4}f(2+\sqrt{2}) \approx 0.02491879 \ll 4,
    \end{equation*}
    and
    \begin{equation*}
      E_2(f) = \frac{f^{(4)}(\tau)}{6} = \frac{4}{(1+\tau)^5} = E_2^{true}(f) \approx 0.02491879,
    \end{equation*}
    which implies $\tau \approx 1.761256$.
    \end{enumerate}
\item
  % 问题4
  \textbf{Solution}\par
  \begin{enumerate}
    \item %(a)
    Notice that $\forall i = 1,2,\dots,n$, there exist
    \begin{align*}
      p(x_i) &= \sum_{m=1}^{n}[h_m(x_i)f_m + q_m(x_i)f'_m] = (a_i+b_ix_i)f_i + (c_i+d_ix_i)f_i' = f_i\\
      p'(x_i) &= \sum_{m=1}^n[h_m'(x_i)f_m+q_m'(x_i)f_m'] = [b_i+(a_i+b_ix_i)2l_i'(x_i)]f_i + [d_i+(c_i+d_ix_i)2l_i'(x_i)]f_i' = f_i'
    \end{align*}
    Set
    \begin{equation*}
      \begin{cases}
        a_i + b_ix_i &= 1\\
        c_i + d_ix_i &= 0\\
        b_i + (a_i+b_ix_i)2l_i'(x_i) &= 0\\
        d_i + (c_i+d_ix_i)2l_i'(x_i) &= 1
      \end{cases}
      \forall i = 1,2,\dots,n.
    \end{equation*}
    Then we get
    \begin{equation*}
      a_m = 1 + 2x_ml'_m(x_m), b_m = -2l_m'(x_m), c_m  = -x_m,d_m = 1, \forall m = 1,2,\dots,n.
    \end{equation*}
  \item %(b)
    By the definition of Hermite interpolation, $\forall f \in \mathbb{P}_{2n-1}$, p(f) = f. Hence we have
    \begin{align*}
      I_n(f) = I_n(p(f)) &= \int\rho \sum_{k=1}^n(h_kf_k+q_kf'_k)\\
      &= \sum_{k=1}^n[f(x_k)\int\rho h_k+f'(x_k)\int \rho q_k]
      & = \sum_{k=1}^n[\omega_kf(x_k)+\mu_kf'(x_k)],
    \end{align*}
    where
    \begin{align*}
      \omega_k &= \int\rho(t)h_k(t)dt = \int\rho_t(1+2x_kl_k'(x_k)-2l_k'(x_k)t)l_k^2(t)dt,\\
      \mu_k &= \int\rho(t)q_k(t)dt = \int\rho(t)(t-x_k)l_k^2(t)dt.
    \end{align*}
  \item %(c)
    Let
    \begin{align*}
      \mu_k &= \int\rho(t)(t-x_k)l_k^2(t)dt = 0\\
      & \implies \int\rho(t)v_n(t)l_k(t) = 0,\forall k = 1,2,\dots,n.
    \end{align*}
    Use the fact that $span\{l_1(t),l_2(t),\dots,l_n(t)\} = span\{1,t,\dots,t^{n-1}\}$, and we have
    \begin{equation*}
      \int \rho(t)v_n(t)p(t)dt = 0, \forall p\in \mathbb{P}_{n-1}.
    \end{equation*}    
    \end{enumerate}

\item
  % 问题5
  \begin{proof}
    For $D^2u(\bar{x})$,
    \begin{align*}
      &\frac{u(\bar{x}-h)-2u(\bar{x})+u(\bar{x}+h)}{h^2} - u''(\bar{x}) = \\
      &\frac{u(\bar{x})-hu'(\bar{x})+\frac{1}{2}h^2u''(\bar{x})-\frac{1}{6}h^3u'''(\bar{x})+\frac{1}{24}h^4u^{(4)}(\bar{x})+O(h^5)}{h^2} +\\
      & \frac{u(\bar{x})+hu'(\bar{x})+\frac{1}{2}h^2u''(\bar{x})+\frac{1}{6}h^3u'''(\bar{x})+\frac{1}{24}h^4u^{(4)}(\bar{x})+O(h^5)}{h^2} - \frac{2u(\bar{x})}{h^2} - u''(\bar{x})\\
      &= \frac{1}{12}h^2u^{(4)}(\bar{x}) + O(h^3),
    \end{align*}
    which implies that $D^2u(\bar{x})$ is second-order accurate.\par
    Next, denote the actual input value with random erors $\epsilon$ by $\widetilde{u}$, and the exact value by $u$. Then
    \begin{equation*}
      \vert u''(\bar{x}) - D^2u(\bar{x})\vert \leq \vert u''(\bar{x})-D_T^2(u(\bar{x})) \vert + \vert D_T^2u(\bar{x})-D^2u(\bar{x})\vert.
    \end{equation*}
    The first term on the right-hand side can be derived by above result:
    \begin{equation*}
       \vert u''(\bar{x}) - D_T^2u(\bar{x})\vert = \vert u''(\bar{x}) - \frac{u(\bar{x}-h)-2u(\bar{x})+u(\bar{x}+h)}{h^2} = \frac{h^2}{12}\vert u^{(4)}(\xi)\vert,
    \end{equation*}
    where $\xi \in [\bar{x}-h,\bar{x}+h]$. \par
    For the second item,
    \begin{align*}
      \vert D_T^2u(\bar{x}) - D^2u(\bar{x})\vert &= \vert \frac{u(\bar{x}-h)-2u(\bar{x})+u(\bar{x}+h)}{h^2} - \frac{\widetilde{u}(\bar{x}-h)-2\widetilde{u}(\bar{x})+\widetilde{u}(\bar{x}+h)}{h^2}\\
      & \leq \frac{1}{h^2}(\vert u(\bar{x}-h)- \widetilde{u}(\bar{x}-h) \vert + 2\vert u(\bar{x}) - \widetilde{u}(\bar{x})\vert + \vert u(\bar{x}+h) - \widetilde{u}(\bar{x}+h)\vert )\\
      \leq \frac{4E}{h^2}.
    \end{align*}
    Therefore we have
    \begin{equation*}
      \vert u''(\bar{x}) - D^2u(\bar{x})\vert \leq \frac{h^2}{12}\vert u^{(4)}(\xi)\vert + \frac{4E}{h^2},
    \end{equation*}
    where $\xi \in[\bar{x}-h,\bar{x}+h]$.\par
    Notice that
    \begin{equation*}
      \frac{h^2}{12}\vert u^{(4)}(\xi)\vert + \frac{4E}{h^2} \geq 2\sqrt{\frac{h^2}{12}\vert u^{(4)}(\xi)\vert \cdot \frac{4E}{h^2}} = 2\sqrt{\frac{E\vert u^{(4)}(\xi)\vert}{3}},
    \end{equation*}
    where $\frac{h^2}{12}\vert = \frac{4E}{h^2}$ yields $h=(\frac{48E}{\vert u^{(4)}(\xi)\vert})^{\frac{1}{4}}$.\par
    Lastly, we design a fourth-order accurate formula based on a symmetric stencil as follows:
    
    \end{proof}
  
\end{enumerate}



\end{document}
